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1 Introduction 



Current and planned experiments in high energy physics can probe physics in 
processes with polarized beams and many tagged particles in the final state. 
The combinatorial explosion of the number of Feynman diagrams contribut- 
ing to scattering amplitudes for many external particles calls for the devel- 
opment of more compact representations that translate well to efficient and 
reliable numerical code. In gauge theories, the contributions from individ- 
ual Feynman diagrams are gauge dependent. Strong numerical cancellations 
in a redundant representation built from individual Feynman diagrams lead 
to a loss of numerical precision, stressing further the need for eliminating 
redundancies. 

Due to the large number of processes that have to be studied in order 
to unleash the potential of modern experiments, the construction of nearly 
optimal representations must be possible algorithmically on a computer and 
should not require human ingenuity for each new application. 

O'Mega JI], H, ^ is a compiler for tree- level scattering amplitudes that 
satisfies these requirements. O'Mega is independent of the target language 
and can therefore create code in any programming language for which a 
simple output module has been written. To support a physics model, O'Mega 
requires as input only the Feynman rules and the relations among coupling 
constants. 

Similar to the earlier numerical approaches and 0, O'Mega reduces 
the growth in calculational effort from a factorial of the number of particles to 
an exponential. The symbolic nature of O'Mega, however, increases its flex- 
ibility. Indeed, O'Mega can emulate both and and produces code that 
is empirically at least twice as fast. The detailed description of all algorithms 
is contained in the extensively commented source code of O'Mega |Q. 

In this note, we sketch the architecture of O'Mega and describe the usage 
of the flrst version. The building blocks of the representation of scattering 
amplitudes generated by O'Mega are described in section |^ and directed 
acyclical graphs are introduced in section 0. The algorithm for constructing 
the directed acyclical graph is presented in section ^ and its implementation 
is described in section |^. We conclude with a few results and examples in 
section Practical information is presented in the appendices: installation 
of the O'Mega software in appendix ^ running of the O'Mega compiler in 
appendix |B| and using O'Mega's output in appendix 0. Finally, appendix 
briefly discusses mechanisms for extending O'Mega. 
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2 One Particle Off Shell Wave Functions 



One Particle Off-Shell Wave Functions (IPOWs) are obtained from con- 
nected Greensfunctions by applying the LSZ reduction formula to all but 
one external line while the remaining line is kept off the mass shell 

,0(g^);out|$(a;)|0(pi), . . . ,0(p„);in) . (1) 

Depending on the context, the off shell line will either be understood as am- 
putated or not. For example, 0(^2); out|<l>(a;) in) in unflavored 
scalar (p^-theoij is given at tree level by 



X X X X 




Pi qi Pi qi Pi qi Pi qi 

The number of distinct momenta that can be formed from n external 
momenta is P{n) = 2"~^ — 1. Therefore, the number of tree IPOWs grows 
exponentially with the number of external particles and not with a factorial, 
as the number of Feynman diagrams, e.g. F{n) = {2n — 5)!! = {2n — 5) ■ 
... 5 ■ 3 ■ 1 in unflavored 0^-theory. 

At tree-level, the set of all IPOWs for a given set of external momenta 
can be constructed recursively 

X X 

where the sum extends over all partitions of the set of n momenta. This 
recursion will terminate at the external wave functions. 

For all quantum field theories, there are — well defined, but not unique — 
sets of Keystones K such that the sum of tree Feynman diagrams for a 
given process can be expressed as a sparse sum of products of IPOWs without 
double counting. In a theory with only cubic couplings this is expressed as 

F{n) P(n) 

^=E^^= E Klf^}SPk,PUPrn)WfMWfXpi)WUPrn), (4) 
i=l k,l,m=l 

with obvious generalizations. The non-trivial problem is to avoide the double 
counting of diagrams like 
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where the circle denotes the keystone. The problem has been solved ex- 
plicitely for general theories with vertices of arbitrary degrees [0. The solu- 
tion is inspired by arguments [Q based on the equations of motion (EOM) 
of the theory in the presence of sources. The iterative solution of the EOM 
leads to the construcion of the IPOWs and the constraints imposed on the 
IPOWs by the EOM suggest the correct set of partitions {{pk,Pi,Pm)} in 
equation (^). 

The maximally symmetric solution selects among equivalent diagrams the 
keystone closest to the center of a diagram. This corresponds to the numerical 
expressions of 0]. The absence of double counting can be demonstrated by 
counting the number F{dma,x, n) of unflavored Feynman tree diagrams with n 
external legs and vertices of maximum degree c/max in to different ways: once 
directly and then as a sum over keystones. The number F^dmax, ^d,n) of 
unflavored Feynman tree diagrams for one keystone N^^n = {^i, n2, . . . , n^}, 
with n = rii + 77-2 + ■ ■ ■ + n^, is given by the product of the number of subtrees 
and symmetry factors 

AT \ n! F{d^^^,ni + 1) 
F{d^^^, Nd,n) = , r I I j 5a 

where |iS(A^)| is the size of the symmetric group of A^, a{n,2n) = 2 and 
a{n, m) = 1 otherwise. Indeed, it can be verified that the sum over all 
keystones reproduces the number 

F(c^max,n) = 5^ J2 ^Kax,iV) (5b) 

d=3 N={ni,n2,... ,nii} 

ni+n2A \-na=n 

l<ni<n2<---<njj<[n/2j 

of all unflavored Feynman tree diagrams. 

A second consistent prescription for the construction of keystones is max- 
imally asymmetric and selects the keystone adjacent to a chosen external 
line. This prescription reproduces the approach in where the tree-level 
Schwinger-Dyson equations are used as a special case of the EOM. 

Recursive algorithms for gauge theory amplitudes have been pioneered 
in 0. The use of IPOWs as basic building blocks for the calculation of 
scattering amplitudes in tree approximation has been advocated in and a 
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heuristic procedure, without reference to keystones, for minimizing the num- 
ber of arithmetical operations has been suggested. This approach is used by 
MADGRAPH ^ for fully automated calculations. The heuristic optimiza- 
tions are quite efficient for 2 —>■ 4 processes, but the number of operations 
remains bounded from below by the number of Feynman diagrams. 



2.1 Ward Identities 

A particularly convenient property of the IPOWs in gauge theories is that, 
even for vector particles, the IPOWs are 'almost' physical objects and satisfy 
simple Ward Identities 

d 

(out|A^(x)|in) =0 (6a) 



Qrf, ^ I ^ ' I ' amp. 



for unbroken gauge theories and 
d 



dXf, ^^^^l^^^^^l^^'^amp. = -^w (out|0H^(a;)|in)^^p (6b) 

for spontaneously broken gauge theories in i?g-gauge for all physical external 
states |m) and \out). Thus the identities (|^) can serve as powerful numer- 
ical checks both for the consistency of a set of Feynman rules and for the 
numerical stability of the generated code. The code for matrix elements can 
optionally be instrumented by O'Mega with numerical checks of these Ward 
identities for intermediate lines. 



3 Directed Acyclical Graphs 

The algebraic expression for the tree-level scattering amplitude in terms of 
Feynman diagrams is itself a tree. The much slower growth of the set of 
IPOWs compared to the set of Feynman diagrams shows that this repre- 
sentation is extremely redundant. In this case. Directed Acyclical Graphs 
(DAGs) provide a more efficient representation, as illustrated by a trivial 
example 

where one multiplication is saved. The replacement of expression trees by 
equivalent DAGs is part of the repertoire of optimizing compilers, known 
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as common subexpression elimination. Unfortunately, this approach fails in 
practice for all interesting expressions appearing in quantum field theory, 
because of the combinatorial growth of space and time required to find an 
almost optimal factorization. 

However, the recursive definition in equation (|^) allows to construct the 
DAG of the IPOWs in equation (^) directly [0, without having to construct 
and factorize the Feynman diagrams explicitely. 

As mentioned above, there is more than one consistent prescription for 
constructing the set of keystones |jl|] . The symbolic expressions constructed 
by O'Mega contain the symbolic equivalents of the numerical expressions 
computed by (maximally symmetric keystones) and (maximally asym- 
metric keystones) as special cases. 

4 Algorithm 

By virtue of their recursive construction in Eqs. (^, tree-level IPOWs form 
a DAG and the problem is to find the smallest DAG that corresponds to 
a given tree, (i.e. a given sum of Feynman diagrams). O'Mega's algorithm 
proceeds in four steps 

Grow: starting from the external particles, build the tower of all 
IPOWs up to a given height (the height is less than the number 
of external lines for asymmetric keystones and less than half of 
that for symmetric keystones) and translate it to the equivalent 
DAG D. 

Select: from D, determine all possible flavored keystones for 
the process under consideration and the IPOWs appearing in 
them. 

Harvest: construct a sub- DAG D* C D consisting only of nodes 
that contribute to the IPOWs appearing in the flavored key- 
stones. 

Calculate: multiply the IPOWs as specified by the keystones 
and sum the keystones. 

By construction, the resulting expression contains no more redundancies and 
can be translated to a numerical expression. In general, asymmetric key- 
stones create an expression that is smaller by a few percent than the result 
from symmetric keystones, but it is not yet clear which approach produces 
the numerically more robust results. 
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The details of this algorithm as implemented in O'Mega are described in 
the source code 0]. The persistent data structures used for the determi- 
nation of D* are very efficient so that the generation of, e. g. Fortran code for 
amplitudes in the Standard Model is always much faster than the subsequent 
compilation. 



5 Implementation 

The O'Mega compiler is implemented in O'Caml [TD|, a functional program- 
ming language of the ML family with a very efficient, portable and freely 
available implementation, that can be bootstrapped on all modern comput- 



ers in a few minutes. The library modules built on experience from |TT], |T2 

The powerful module system of O'Caml allows an efficient and concise 
implementation of the DAGs for a specific physics model as a functor appli- 
cation [|l|. This functor maps from the category of trees to the category of 
DAGs and is applied to the set of trees defined by the Feynman rules of any 
model under consideration. 

The module system of O'Caml has been used to make the combinato- 
rial core of O'Mega demonstrably independent from the specifics of both the 
physics model and the target language as shown in Figure A For- 
tran90/95 backend has been realized first, backends for C++ and Java will 
follow. The complete electroweak Standard Model has been implemented 
together with anomalous gauge boson couplings. The implementation of in- 
terfering color amplitudes is currently being completed. 

Many extensions of the Standard Model, most prominently the Minimal 
Supersymmetric Standard Model (MSSM), contain Majorana fermions. In 
this case, fermion lines have no canonical orientation and the determination 
of the relative signs of interfering amplitudes is not trivial. However, the 
Feynman rules for Majorana fermions and fermion number violating interac- 
tions proposed in [|1^ have been implemented in O'Mega in analogy to the 
naive Feynman rules for Dirac fermions and both methods are available. Nu- 
merical comparisons of amplitudes for Dirac fermions calculated both ways 
show agreement at a small multiple of the machine precision. Thus, all in- 
gredients for the MSSM are available in O'Mega and the implementation of 
the complete MSSM Lagrangian is currently under way. Non-minimal gauge 
models, including left-right symmetric models, can be implemented easily. 

As mentioned above, the compilers for the target programming language 
are the slowest step in the generation of executable code. On the other 
hand, the execution speed of the code is limited by non-trivial vertex eval- 
uations for vectors and spinors, which need 0(10) complex multiplications. 
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Figure 1: Module dependencies in O'Mega. The diamond shaped nodes 
denote abstract signatures defining functor domains and co-domains. The 
rectangular boxes denote modules and functors, while oval boxes stand for 
example applications. 
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Table 1: Radiative corrections to four fermion production: comparison of 
the computational complexity of scattering amplitudes obtained from Feyn- 
man diagrams and from O'Mega. (The counts correspond to the full Stan- 
dard Model — sans light fermion Yukawa couplings — in unitarity gauge with 
quartic couplings emulated by cubic couplings of non-propagating auxiliary 
fields.) 



process 

e^e^ —>■ e^Uedubb 
e^e^ — > e'^h'edubb'y 
e'^e~ — >• e'^Uedubb'y-)' 



Diagrams 
^ vertices 
472 2832 
4956 34692 
58340 466720 



O'Mega 
#prop. vertices 

49 232 
108 722 
226 2212 



Table 2: Radiative corrections to six fermion production: comparison of 
the computational complexity of scattering amplitudes obtained from Feyn- 
man diagrams and from O'Mega. (The counts correspond to the full Stan- 
dard Model — sans light fermion Yukawa couplings — in unitarity gauge with 
quartic couplings emulated by cubic couphngs of non-propagating auxiliary 
fields.) 



Therefore, an O'Mega Virtual Machine can challenge native code and avoid 
compilations. 



6 Results 
6.1 Examples 

Tables |I] and |^ show the reduction in computational complexity for some 
important processes at a e"'"e~-linear collider including radiative corrections. 
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Using the asymmetric keystones can reduce the number of vertices by some 10 
to 20 percent relativ to the quoted numbers for symmetric keystones. 



6.2 Comparisons 

HELAC's [^] diagnostics report more vertices than O'Mega for identical am- 
phtudes. This ranges from comparable numbers for Standard Model pro- 
cesses with many different flavors to an increase by 50 percent for processes 
with many identical flavors. Empirically, O'Mega's straight line code is twice 
as fast as HELAC's DO-loops for identical optimizing Fortran95 compilers 
(not counting HELAC's initialization phase). Together this results in an 
improved performance by a factor of two to three. 

The numerical efficiency of O'Mega's Fortran95 runtime library is empir- 
ically identical to HELAS j^. Therefore, O'Mega's performance can directly 
be compared to MADGRAPH's by comparing the number of vertices. For 
2 — > 5-processes in the Standard Model, O'Mega's advantage in performance 
is about a factor of two and grows from there. 

The results have been compared with MADGRAPH |^ for many Stan- 
dard Model processes and numerical agreement at the level of 10~^^ has been 
found with double precision floating point arithmetic. 



6.3 Applications 

O'Mega generated amplitudes are used in the omnipurpose event generator 
generator WHIZARD The first complete experimental study of vector 
boson scattering in six fermion production for linear collider physics [|T5| has 
been facilitated by O'Mega and WHIZARD. 
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A Installing O'Mega 
A.l Sources 



O'Mega is Free Software and the sources can be obtained from [http : / / 



heplix ■ ikp . physik . tu-darmstadt . de/~ohl/ omega/] or from [ftp : //heplix 



ikp ■ physik . tu-darmstadt . de/pub/ ohl/omeg4 The command 



ohl(§thopad: ~mc$ zcat omega-current .tar .gz I tar xf - 

will unpack the sources to the directory |omega| . The subdirectories of [omega 



are 



bin contains executable instances of O'Mega: f90_SM.bin (f 90_ 



SM . opt| if the sytem is supported by O'Caml's native code com- 



piler), |f90_qED.binl , etc. 

doc contains I^TgX sources of user documentation. 

examples contains currently no supported examples. 

lib contains library support for targets (Fortran90/95 modules, 
etc.). 

src contains the unabridged and uncensored sources of O'Mega, 
including comments. 

tests contains a battery of regression tests. Most tests require 
Madgraph [§]. 

web contains the 'woven' sources, i. e. a pretty printed version 
of the source including BTgX documentation. Weaving the 
sources requires programs, |ocamlweb| and |noweb| . A complete 



PostScript file is available from the same place as the O'Mega 
sources. (It is not required for the end user to read this.) 
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A. 2 Prerequisites 

A. 2.1 Objective Caml (a. k. a. O'Caml) 

You need version 3.00 or liigtierQ. You can get it from [http : //pauillaic" 



inria . f r/ocaiiil7| . There are precompiled binaries for some popular systems 
and complete sources. Building from source is straightforward (just follow 
the instructions in the file [INSTALL in the toplevel directory, the defaults 



are almost always sufficient) and takes (9(10) minutes on a modern desk- 
top system. If available for your system (cf. the file [README] in the toplevel 
directory), you should build the native code compiler. 

A. 2. 2 GNU make 

This should be available for any system of practical importance and it makes 
no sense to waste physicist's time on supporting all incompatible flavors of 



nake in existence. GNU make is the default on Linux systems and is often 



available as |gmake| on commercial Unices. 



A. 2. 3 Fortran90/95 Compiler 

Not required for compiling or running O'Mega, but Fortran90/95 is currently 
the only fully supported target. 

A. 3 Configuration 

Before the next step, O'Caml must have been installed. Configuration is 
performed automatically by testing some system features with the command 

$ ./configure 

See 

$ ./configure — help 



for additional options. NB: The use of the options | — enable-gui| and 



-enable-unsupported] is strongly discouraged. The resulting programs 
require additional prerequisites and even if you can get them to compile, the 
results are unpredictable and we will not answer any questions about them. 



NB: |conf igure] keeps it's state in |conf ig. cach4 If you want to reconfigure 



after adding new libraries to your system, you should remove [config. cache 



before running |conf igure 



^O'Mega can probably be made to work with versions 2. Ox. Given the simphcity of 
building O'Caml version 3.00, this is not worth the effort. 
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A. 4 Compilation 

The command 



$ make bin 

will build the byte code executables. For each pairing of physics model and 
target language, there will be one executable. 

$ make opt 

will build the native code executables if the sytem is supported by O'Caml's 
native code compiler and it is installed. The command 

$ make f95 

will build the Fortran90/95 hbrary and requires, obviously, a Fortran90/95 
compiler. 



B Running O'Mega 

O'Mega is a simple application that takes parameters from the commandline 
and writes results to the standard output device^ (diagnostics go to the 
standard error device). E.g., the UNIX commandline 

$ ./bin/f90_SM.opt e+ e- e+ nue ubar d > cc20_amplitude . f 95 

will cause O'Mega to write a Fortran95 module containing the Standard 



Model tree level scattering amplitude for e^e e^Ueud to the file |cc20_ 



amplitude . f 95| . Particles can be combined with colons. E.g., 

$ ./bin/f90_SM.opt ubar : u : dbar : d ubar : u : dbar : d e+:mu+ e-:mu- > dy.f95 

will cause O'Mega to write a Fortran95 module containing the Standard 
Model tree level parton scattering amplitudes for all Drell-Yan processes to 



the file |dy.f95 



A synopsis of the available options, in particular the particle names, can 
be requested by giving an illegal option, e. g.: 



^In the future, other targets than Fortran90/95 might require more than one output 
file (e. g. source files and header files for C/C++). In this case the filenames will be specified 
by commandline parameters. 



14 



$ ./bin/f90_SM.opt -? 
./bin/f90_SM.opt: unknown option '-?'. 

usage: . /bin/f 90_SM. opt [options] [e- I nue I u I d I e+ I nuebar I ubar I dbar\ 
I mu- I numu I c I s I mu+ I numubar I cbar I sbar I tau- I nutau I t I b\ 
I tau+ 1 nutaubar I tbar I bbar I A I Z I W+ 1 W- 1 g I H I phi+ I phi- I phiO] 
-target : function function name 

-target: 90 don't use Fortran95 features that are not in Fortran90 
-target: kind real and complex kind (default: omega_prec) 
-target : width approx. line length 
-target : module module name 
-target: use use module 

-target :whizard include WHIZARD interface 

-model : constant_width use constant width (also in t-channel) 

-model : fudged_width use fudge factor for charge particle width 

-model : custom_width use custom width 

-model : cancel_widths use vanishing width 

-warning: check arguments and print warning on error 

-error: check arguments and terminate on error 

-warning: a check # of input arguments and print warning on error 
-error: a check # of input arguments and terminate on error 
-warning :h check input helicities and print warning on error 
-error :h check input helicities and terminate on error 
-warning :m check input momenta and print warning on error 
-error :m check input momenta and terminate on error 
-warning :g check internal Ward identities and print warning on error 
-error :g check internal Ward identities and terminate on error 
-forest ??? 

-revision print revision control information 

-quiet don't print a summary 

-summary print only a summary 

-params print the model parameters 

-poles print the Monte Carlo poles 

-dag print minimal DAG 

-full_dag print complete DAG 

-file read commands from file 



B.l General Options 

-warning : include code that checks the supplied arguments and 
prints a warning in case of an error. 
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-warning : a check the number of input arguments (momenta and 
spins) and print a warning in case of an error. 

-warning :h check the values of the input hehcities and print a 
warning in case of an error. 

-warning :m check the values of the input momenta and print a 
warning in case of an error. 

-warning :g check internal Ward identities and print a warning 
in case of an error (not supported yet!). 

-error: like -warning: but terminates on error. 

-error: a like -warning: a but terminates on error. 

-error :h like -warning :h but terminates on error. 

-error :m like -warning :m but terminates on error. 

-error :g like -warning :g but terminates on error. 

-revision print revision control information 

-quiet don't print a summary 

-summary print only a summary 

-params print the model parameters 

-poles print the Monte Carlo poles in a format understood by 
the WHIZARD program [|l|]. 

-dag print the reduced DAG in a format understood by the dot 
program. 

-full_dag print the complete DAG in a format understood by 
the dot program. 

-file read commands from file 
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B.2 Model Options 
B.2.1 Standard Model 

-model : constant_width use constant width (also in i-channel) 

-model : f udged_width use fudge factor for charge particle width 

-model : custom_width use custom width 

-model : cancel-widths use vanishing width 

B.3 Target Options 

B. 3.1 Fortran90/95 

-target : function function name 

-target : 90 don't use Fortran95 features that are not in For- 
tran90 

-target: kind real and complex kind (default: omega_prec) 
-target : width approx. line length 
-target : module module name 
-target: use use module 

-target :whizard include WHIZARD interface 

C Using O'Mega's Output 

The structure of the outputfile, the calling convention and the required li- 
braries depends on the target language, of course. 

C. l Fortran90/95 

The Fortran95 module written by O'Mega has the following signature 
module omega_ amplitude 
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C.1.1 Libraries 



The imported Fortran modules are 

omega_kinds defines omega_prec, which can be whatever the For- 
tran compiler supports. NB: the support libraries have not yet 
been tuned to give reliable answers for amplitudes with gauge 
cancellations in single precision. 

omega95 defines the vertices for Dirac spinors in the chiral repre- 
sentation and vectors. 

omega95_bispinors is an alternative that defines the vertices for 
Dirac and Majorana spinors in the chiral representation and 
vectors using the Feynman rules of [pT^ . 



omega_parameters defines the coupling constants 



use omega_kinds 

use omega95 

use omega_parameters 

implicit none 

private 



C.1.2 Summary of Exported Functions 

The functions and subroutines experted by the Fortran95 module are 

• the scattering amplitude in different flavor bases (arrays of PDG codes 
or internal numbering): 

public :: amplitude, ainplitude_f , amplitude_l, amplitude_2 

• the scattering amplitude with heuristics supressing vanishing helicity 
combinations: 

public : : amplitude_nonzero , amplitude_f _nonzero , & 
amplitude_l_nonzero , amplitude_2_nonzero 

• the squared scattering amplitude summed over helicity states 
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public :: spin_sum_sqme , spin_sum_sqme_l , suin_sqme 
public : : spin_sum_sqme_nonzero , spin_suin_sqme_l_nonzero , & 
sum_sqme_nonzero 



• "scattering" a general density matrix 

public : : scatter, scatter_nonzero 

• "scattering" a diagonal density matrix 

public :: scatter_diagonal , scatter_diagonal_nonzero 

• inquiry and maintenance functions 

public : : allocate_zero 

public : : multiplicities, multiplicities_in, multiplicities_out 

public :: number_particles , & 

number_particles_in, number_particles_out 

public :: number_spin_states , & 

number_spin_states_in, number_spin_states_out , & 
spin_states, spin_states_in, spin_states_out 

public :: number_f lavor_states , & 

number_f lavor_states_in, number_f lavor_states_out , & 
f lavor_states , f lavor_states_in, f lavor_states_out 

public : : number_f lavor_zeros , & 

number_f lavor_zeros_in, number_f lavor_zeros_out , & 
f lavor_zeros , f lavor_zeros_in, f lavor_zeros_out 

public : : create, reset, destroy 

C.1.3 Maintenance Functions 



They currently do nothing, but are here for WHIZARD's \T^] convenience 
create is called only once at the very beginning, 
reset is called whenever parameters are changed, 
destroy is called at most once at the very end. 
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subroutine create () 
end subroutine create 

subroutine reset () 
end subroutine reset 
subroutine destroy () 
end subroutine destroy 

Allocate an array of the size used by the heuristic that suppresses vanishing 
helicity combinations 

interface allocate_zero 

module procedure allocate_zero_l , allocate_zero_2 
end interface 

for join numbering of in and out states 

subroutine allocate_zero_l (zero) 

integer, dimension( : , : ) , pointer :: zero 
end subroutine allocate_zero_index 

and for separate numbering of in and out states 

subroutine allocate_zero_2 (zero) 

integer, dimension (: ,:,:,:), pointer : : zero 
end subroutine allocate_zero_index_inout 

C.1.4 Inquiry Functions 

The total number of particles, the number of incoming particles and the 
number of outgoing particles: 

pure function number_particles () result (n) 

integer : : n 
end function number_particles 

pure function number_particles_in () result (n) 

integer : : n 
end function number_particles_in 
pure function number_particles_out () result (n) 

integer : : n 
end function number_particles_out 

The spin states of all particles that can give non-zero results and their num- 
ber. The tables are interpreted as 



20 



s (1 : , i) contains the helicities for each particle for the ith hehc- 
ity combination. 

pure function number. spin_ states () result (n) 

integer : : n 
end function number_spin_states 
pure subroutine spin_states (s) 

integer, dimension( : , : ) , intent (inout) :: s 
end subroutine spin_states 

The spin states of the incoming particles that can give non-zero results and 
their number: 

pure function number_spin_states_in () result (n) 
integer : : n 

end function number_spin_states_in 
pure subroutine spin_states_in (s) 

integer, dimension( : , : ) , intent (inout) :: s 
end subroutine spin_states_in 

The spin states of the outgoing particles that can give non-zero results and 
their number: 

pure function number_spin_states_out () result (n) 

integer : : n 
end function number_spin_states_out 
pure subroutine spin_states_out (s) 

integer, dimension( : , : ) , intent (inout) :: s 
end subroutine spin_states_out 

The flavor combinations of all particles that can give non-zero results and 
their number. The tables are interpreted as 

f (1 : , i) contains the PDG particle code for each particle for the 
ith helicity combination. 

pure function number_f lavor_states () result (n) 

integer : : n 
end function number_f lavor_states 
pure subroutine f lavor_states (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine f lavor_states 
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The flavor combinations of the incoming particles that can give non-zero 
results and their number. 

pure function number_f lavor_states_in () result (n) 

integer : : n 
end function number.f lavor_states_in 
pure subroutine f lavor_states_in (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine f lavor_states_in 

The flavor combinations of the outgoing particles that can give non-zero 
results and their number. 

pure function number_f lavor_states_out () result (n) 

integer : : n 
end function number_f lavor_states_out 
pure subroutine f lavor_states_out (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine f lavor_states_out 

The flavor combinations of all particles that always can give a zero result 
and their number: 

pure function number_f lavor_zeros () result (n) 

integer : : n 
end function number_f lavor_zeros 
pure subroutine flavor_zeros (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine flavor_zeros 

The flavor combinations of the incoming particles that always can give a zero 
result and their number: 

pure function number_f lavor_zeros_in () result (n) 

integer : : n 
end function number_f lavor_zeros_in 
pure subroutine f lavor_zeros_in (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine f lavor_zeros_in 

The flavor combinations of the outgoing particles that always can give a zero 
result and their number: 
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pure function number_f lavor_zeros_out () result (n) 

integer : : n 
end function number_f lavor_zeros_out 
pure subroutine f lavor_zeros_out (f) 

integer, dimension( : , : ) , intent (inout) :: f 
end subroutine f lavor_zeros_out 

The same initial and final state can appear more than once in the tensor 
product and we must avoid double counting. 

pure subroutine multiplicities (a) 

integer, dimension( : ) , intent (inout) :: a 
end subroutine multiplicities 

pure subroutine multiplicities_in (a) 

integer, dimension( : ) , intent (inout) :: a 
end subroutine multiplicities_in 

pure subroutine multiplicities_out (a) 

integer, dimension( : ) , intent (inout) :: a 
end subroutine multiplicities_out 



C.1.5 Amplitude 

The function arguments of of the amplitude are 

k(0:3,l:) are the particle momenta: k(0:3,l) and k(0:3,2) 
are the incoming momenta, k (0:3,3:) arc the outgoing mo- 
menta. All momenta are the physical momenta, i.e. forward 
time-like or light-like. The signs of the incoming momenta are 
flipped internally. Unless asked by a commandline parameter, 
O'Mega will not check the validity of the momenta. 

s (1 : ) are the hehcities in the same order as the momenta, s = ±1 
signify s — ±1/2 for fermions. s — makes no sense for 
fermions and massless vector bosons s = 4 signifies an un- 
physical polarization for vector boson that the users are not 
supposed to use. Unless asked by a commandline parameter, 
O'Mega will not check the validity of the hehcities. 

f (1:) are the PDG particle codes in the same order as the mo- 
menta. 
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pure function amplitude (k, s, f) result (amp) 

real(kind=omega_prec) , dimension (0 :,:) , intent (in) :: k 
integer, dimension( : ) , intent(in) : : s, f 
complex (kind=omega_prec) : : amp 

end function amplitude 

Identical to amplitude (k, s, flavors(: ,f)), where flavors has been 
filled by f lavor_states: 

pure function amplitude_f (k, s, f) result (amp) 

real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 

integer, dimension( : ) , intent(in) :: s 

integer, intent (in) : : f 

complex (kind=omega_prec) : : amp 
end function amplitude_f 

Identical to amplitude (k, spins(:,s), f lavors( : ,f)), where spins has 
been filled by spin_states and flavors has been filled by f lavor_states: 

pure function amplitude_l (k, s, f) result (amp) 

real(kind=omega_prec) , dimension (0 :,:) , intent (in) :: k 
integer, intent (in) : : s, f 
complex (kind=omega_prec) : : amp 

end function amplitude_l 

Similar to amplitude. 1, but with separate incoming and outgoing particles: 

pure function amplitude_2 & 

(k, s_in, f_in, s_out, f_out) result (amp) 
real(kind=omega_prec) , dimension (0: ,:) , intent (in) :: k 
integer, intent (in) : : s_in, f_in, s_out, f_out 
complex (kind=omega_prec) : : amp 

end function amplitude_2 

The following are subroutines and not functions, since Fortran95 restricts ar- 
guments of pure functions to intent (in), but we need to update the counter 
for vanishing amplitudes. 

zerod : , 1 : ) an array containing the number of times a combi- 
nation of spin index and flavor index yielded a vanishing am- 
plitude. After a certain threshold, these combinations will be 
skipped. allocate_zero will allocate the correct size. 

n the current event count 
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pure subroutine amplitude.nonzero (amp, k, s, f, zero, n) 
complex (kind=omega_prec) , intent (out) :: amp 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, dimension( : ) , intent(in) :: s, f 
integer, dimension( : , : ) , intent (inout) :: zero 
integer, intent (in) : : n 

end subroutine amplitude_nonzero 



pure subroutine amplitude_l_nonzero (amp, 
complex (kind=omega_prec) , intent(out) : 
real (kind=omega_prec) , dimension(0 : , : ) , 
integer, intent (in) : : s, f 
integer, dimension( : , : ) , intent (inout) 
integer, intent (in) : : n 

end subroutine amplitude. l_nonzero 



k, s, f, zero, n) 

amp 

intent (in) : : k 



zero 



pure subroutine amplitude_f _nonzero & 
(amp, k, s, f, zero, n) 
complex (kind=omega_prec) , intent (out) : 
real (kind=omega_prec) , dimension (0 : , : ) , 
integer, dimension( : ) , intent(in) :: s 
integer, intent (in) : : f 
integer, dimension(:, 
integer, intent (in) : 

end subroutine amplitude_f _nonzero 



), intent (inout) 
n 



amp 
intent (in) 



zero 



zero (1 : , 1 : , 1 : , 1 : ) an array containing the number of times a 
combination of incoming and outgoing spin indices and flavor 
indices yielded a vanishing amplitude. allocate_zero will al- 
locate the correct size. 

pure subroutine amplitude_2_nonzero & 

(amp, k, s_in, f_in, s_out, f_out, zero, n) 
complex (kind=omega_prec) , intent (out) :: amp 
real(kind=omega_prec) , dimension (0 :,:) , intent (in) :: k 
integer, intent(in) :: s_in, f_in, s_out, f_out 
integer, dimension ( :,:,:,:), intent (inout) : : zero 
integer, intent (in) : : n 

end subroutine amplitude_2_nonzero 
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C.1.6 Summation 



The the sums of squared matrix elements, the optional mask |smask| can be 
used to sum only a subset of helicities or flavors. 

pure function spin_sum_sqme (k, f, smask) result (amp2) 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, dimension( : ) , intent(in) :: f 
logical, dimensionC : ) , intent (in), optional :: smask 
real(kind=omega_prec) : : amp2 

end function spin_sum_sqme 

pure function spin_sum_sqme_l (k, f, smask) result (amp2) 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, intent (in) : : f 

logical, dimension( : ) , intent (in), optional :: smask 
real(kind=omega_prec) : : amp2 
end function spin_sum_sqme_l 

pure function sum_sqme (k, smask, fmask) result (amp2) 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
logical, dimension( : ) , intent(in), optional :: smask, fmask 
real(kind=omega_prec) : : amp2 

end function sum_sqme 

pure subroutine spin_sum_sqme_nonzero (amp2, k, f, zero, n, smask) 
real(kind=omega_prec) , intent (out) :: amp2 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, dimension( : ) , intent(in) :: f 
integer, dimension( : , : ) , intent (inout) :: zero 
integer, intent (in) : : n 

logical, dimension( : ) , intent (in), optional :: smask 
end subroutine spin_sum_sqme_nonzero 

pure subroutine spin_sum_sqme_l_nonzero (amp2, k, f, zero, n, smask) 
real(kind=omega_prec) , intent (out) :: amp2 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, intent (in) : : f 

integer, dimension( : , : ) , intent (inout) :: zero 
integer, intent (in) : : n 

logical, dimension( : ) , intent (in), optional :: smask 
end subroutine spin_sum_sqme_l_nonzero 
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pure subroutine suin_sqme_nonzero (amp2, k, zero, n, smask, fmask) 
real(kind=omega_prec) , intent (out) :: amp2 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) :: k 
integer, dimension( : , : ) , intent (inout) :: zero 
integer, intent (in) : : n 

logical, dimension( : ) , intent(in), optional :: smask, fmask 
end subroutine sum_sqme_masked_nonzero 

C.1.7 Density Matrix Transforms 

There are also utility functions that implement the transformation of density 
matrices directly 



and avoid double counting 

pure subroutine scatter_correlated (k, rho_in, rho_out) 
real(kind=omega_prec) , dimension (0 :,:) , intent(in) :: k 
complex (kind=omega_prec) , dimension ( :,:,:,:), & 

intent (in) : : rho_in 
complex (kind=omega_prec) , dimension ( :,:,:,:), & 
intent (inout) :: rho_out 
end subroutine scatter_correlated 

pure subroutine scatter_correlated_nonzero & 
(k, rho_in, rho_out, zero, n) 
real(kind=omega_prec) , dimension(0 : , : ) , intent (in) :: k 
complex (kind=omega_prec) , dimension ( :,:,:,:), & 

intent (in) :: rho_in 
complex (kind=omega_prec) , dimension ( :,:,:,:), & 

intent (inout) : : rho_out 
integer, dimension( :,:,:,:) , intent(inout) :: zero 
integer, intent (in) : : n 
end subroutine scatter_correlated_nonzero 



P 



p' = TpT^ 



(8) 




(9) 



27 



In no off-diagonal density matrix elements of the initial state are known, the 
computation can be performed more efficiently: 

p'f = Y,Tf^P^T}, = Y,\Tf^?P^ ^ 

i i 

pure subroutine scatter_diagonal (k, rho_in, rho_out) 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) : 
real(kind=omega_prec) , dimension( : , : ) , intent(in) :: 
real(kind=omega_prec) , dimension( : , : ) , intent (inout) 

end subroutine scatter_diagonal 

pure subroutine scatter_diagonal_nonzero & 
(k, rho_in, rho_out, zero, n) 
real(kind=omega_prec) , dimension(0 : , : ) , intent(in) : 
real(kind=omega_prec) , dimension( : , : ) , intent(in) :: 
real(kind=omega_prec) , dimension( : , : ) , intent (inout) 
integer, dimension(: ,:,:,:), intent (inout) : : zero 
integer, intent (in) : : n 

end subroutine scatter_diagonal_nonzero 

Finis. 

end module omega_amplitude 

NB: the name of the module can be changed by a commandline parameter 
and Fortran95 features like pure can be disabled as well. 

C.2 FORTRAN77 

The preparation of a FORTRAN?? target is straightforward, but tedious and 
will only be considered if there is sufficient demand and support. 

C.3 HELAS 

This target for the HELAS library is incomplete and no longer maintained. 
It was used as an early benchmark for the Fortran90/95 library. No vector 
boson selfcouplings are supported. 



: k 

rho_in 
: : rho_out 



: k 

rho_in 
: : rho_out 
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C.4 C, C++ & Java 



These targets does not exist yet and we solicit suggestions from C++ and 
Java experts on useful calling conventions and suppport libraries that blend 
well with the HEP environments based on these languages. At least one of 
the authors believes that Java would be a better choice, but the political 
momentum behind C++ might cause an early support for C++ anyway. 



D Extending O'Mega 

D.l Adding A New Physics Model 

Currently, this still requires to write O'Caml code. This is not as hard 
as it might sound, because an inspection of [bin/models .mil shows that all 
that is required are some tables of Feynman rules that can easily be written 
by copying and modifyng an existing example, after consulting with |src/ 



Icouplings .mlij or the corresponding chapter in the woven source. In fact. 



having the full power of O'Caml at one's disposal is very helpful for avoiding 
needless repetition. 

Nevertheless, in the near future, there will be some special models that 
can read model specifications from external files. The first one of its kind 
will read CompHEP |]16[ model files. Later there will be a native O'Mega 
model file format, but it will probably go through some iterations. 



D.2 Adding A New Target Language 

This will always require to write O'Caml code, which is again not too hard. 
In addition a library for vertices will be required, unless the target performs 
complete inlining. NB: an early experiment with inlining Fortran proved to 
be an almost complete failure on Linux/Intel PCs. The inlined code was 
huge, absolutely unreadable and only marginally faster. The bulk of the 
computational cost is always in the vertex evaluations and function calls 
create in comparison negligible costs. This observation is system dependent, 
of course, and inlining might be beneficial for other architectures with better 
floating point performance, after all. 
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